Fluctuations of time averages for Langevin dynamics in a binding force field 
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We derive a simple formula for the fluctuations of the time average x(t) around the thermal mean 
(a:) eq for overdamped Brownian motion in a binding potential U(x). Using a backward Fokker- 
Planck equation, introduced by Szabo, Schulten, and Schulten in the context of reaction kinetics, 
we show that for ergodic processes these finite measurement time fluctuations are determined by 
the Boltzmann measure. For the widely applicable logarithmic potential, ergodicity is broken. 
We quantify the large non-ergodic fluctuations and show how they are related to a super-aging 
correlation function. 
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Current technology permits tracking of trajectories of 
individual molecules with exquisite precision. The mo- 
tion of a Brownian particle in a binding potential field 
U{x) is used to model many such physical, biological and 
chemical processes. From statistical mechanics, we know 
that if the process is ergodic, and if the measurement 
time t — > oo, then the time average x{t) = J ' x(t')dt' /t 
is equal to the corresponding ensemble average (a;) cq . In 
experiment the measurement time might be long, but it is 
always finite. Hence it is natural to ask what the fluctua- 
tions of x are. Such an analysis sheds light on deviations 
from the thermal equilibrium average due to finite time 
measurement, a general theme which has attracted much 
interest in the context of fluctuation theorems [l| . The 
Boltzmann measure, due to ergodicity, yields equilibrium 
properties of thermal systems. Surprisingly, we find that 
for Langevin dynamics, the Boltzmann measure also de- 
termines the deviations from ergodicity. 

As we will show, for binding fields U(x) where the 
Fokker-Planck (FP) operator exhibits a discrete eigen- 
spectrum, the fluctuations of the time average x become 
small as time increases, as expected from ordinary er- 
godic statistical mechanics. For this type of field, er- 
godicity is related to the work of Szabo, Schulten, and 
Schulten [2| on the seemingly unrelated problem of reac- 
tion kinetics (see details below) . A more interesting case 
is that of a logarithmic binding field |3[ U (x) ~ Uq m(|a;|) 
when \x\ — > oo, since for such a potential the fluctu- 
ations of x are not small even in the long time limit. 
Here the Boltzmann measure exhibits power law tails, 
P cq (x) oc \x\~ Uo ^ kBT \ Starting at the origin, the parti- 
cle during its evolution tends to sample larger and larger 
values of |x| as illustrated in Fig. [T] Large fluctuations 
in the amplitude of x(t) cause the time average of this 
special process to remain random even in the long time 
limit. In what follows, we calculate the magnitude of 
these fluctuations and show how they are related to a 
super-aging correlation function. Importantly, such log- 
arithmic potentials model many physical systems, rang- 
ing from optical lattices Q , charges in vicinity of a long 
charged polymer Q, DNA dynamics @, membrane in- 
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FIG. 1: (color online) The trajectory of a Brownian particle 
in a logarithmic potential exhibits large amplitude fluctua- 
tions. As a consequence, the time average of the process x 
(red curve) does not converge to a fixed value even though 
Boltzmann equilibrium ensemble average (a;) C q is zero. Here 
U{x) — ln(l + x 2 )/2, ksT — 1/2 and the diffusion constant 
D = l. 



duced forces Q, a nano-particle in a trap Q, to long 
ranged interacting models [9(. At the end of this Let- 
ter we discuss the connection between our theory and a 
recent experiment [Icj . 

Model and observable. Brownian dynamics in a force 
field f(x) = — c\U(x)/dx obeys the equation 
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(i) 



Here 7 is the friction constant, w{t) is Gaussian 
white noise obeying the fluctuation dissipation relation 
{V(t)v(t')) = 2DS(t - t'), and D = k B T/-y according 
to the Einstein relation. From the trajectory x(t) we 
construct the time average x(t) = J* x(t')dt'/t. For a 
binding potential, in the long time limit, x obeys the 
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equilibrium Boltzmann distribution: 
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where Z is the normalizing partition function which is 
assumed to be finite. We consider symmetric potentials 
U(x) — U{—x) and then the ensemble average in equi- 
librium (x) eq = J x oo xP cq (x)dx = 0. If the process is 
crgodic then in the long time limit x — > (x) cq = 0. If 
lim t _j. 00 (x 2 (t)) ^ the process is non-ergodic, where (• • •) 
stands for an ensemble mean. In the second part of our 
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work we show that not all binding potentials satisfy the 
ergodic hypothesis. 

Szabo-Schulten-Schulten equation yields the fluctua- 
tions of the time average. The variance of the time aver- 
age is given by 
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(t)> = ~ j*dt 2 y d*i 



[x(t 2 )x(h)} 
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where (x(t 2 )x(ti)) is the correlation function. For the 
Markovian process under investigation, and for a particle 
starting at the origin at time t = we have [ill ] 



1—2 



(*)> = 



t 2 



dto 



dti 



OO pOO 



oo J — oo 



x 2 x\P (x2,t 2 \xi,ti) P (a;i,£i|0, 0) dx\dx 2 
I 



(4) 



where P{x 2l t 2 \xi 1 t\) is the conditional probability den- 
sity to find the particle on x 2 at time t 2 once it is located 
at x\ at time t\. In the limit of long times, the major 
contribution to the integration over t\ comes from long 
times; hence one replaces P(xi,ii|0,0) with P eq (xi). To 
proceed, it is useful to define 



E{xi,r)dT 



(5) 



where E(x±, r) is the averaged position of a particle at a 
time r after it starts at x\. Two cases are of interest; the 
first is when £(x\) is finite, the other when it diverges. We 
shall start with the former case which is clearly relevant 
to potential fields where the FP eigenspectrum has 
a finite energy gap to the ground state, since then the 
relaxation of E(xi,t) is exponential. From Eq. ((3]) it 
follows that in the long time limit 
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As is well known the backward FP equation [TiJ 
tf FP P(x 2 ,r\x u 0) = JLP(x 2 ,t\ Xi ,0), 



L J FP = D 



d 2 , fjxi) d 
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governs the dynamics where Lp P is the adjoint FP op- 
erator and P(x 2 , 0\xi, 0) = S(x 2 — xi). By definition 
E(xi,t) = j_ x 2 P(x2,T\xi,0)dx 2 which implies 
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L fy> E{x 1 ,t) = —E{x u t) 
with E(xi,(S) — x\. Using Eq. J3]), we find 
L\ v i(x\) = -xx 



(8) 



(9) 



with ^(0) = 0. Eq. (|9|) was obtained previously in [2| in 
the context of reaction kinetics. Eqs. (J3HH]) are so general 
that they could be extended to arbitrary Markovian pro- 
cesses. The latter equations thus serve as a starting point 
for the investigation of fluctuations of time averages for 
a wide class of systems. 

Fluctuations of time averages determined from Boltz- 
mann statistics. Eq. (|9]) is easy to solve, and upon using 
Eq. we find the general formula 
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(10) 

As is well known, Boltzmann statistics can be used to de- 
termine the time average of ergodic processes: x — > (x) 
in the long time limit. Eq. (fT0|) shows that also the 
finite time fluctuations of x are determined by the Boltz- 
mann distribution. Surprisingly, Eq. (|10[) shows that 
the difficult task of finding the entire eigenspectrum 
of the FP operator is not required. Eq. (jTU)) is eas- 
ily generalized to dimensions greater than one, and to 
non-thermal processes whose equilibrium density is non- 
Boltzmannian. As expected from ergodicity, the magni- 
tude of the fluctuations decays to zero with time, pro- 
vided that the integrals in Eq. (fTU|) converge. For ex- 
ample, for the harmonic potential U(x) — muj 2 x 2 /2 we 
get (x 2 ) - 2(k B T) 2 /{D(muj 2 ) 2 t}. An interesting case 
where the integrals diverge is the logarithmic potential 
U(x) - U ln{\x\) for |af| ->■ oo and U /(k B T) < 5. This 
leads to a non-ergodic behavior which we now investigate. 

Logarithmic potential. We will first find the two-point 
correlation function (x(t 2 )x(ti)) for a general logarithmic 
potential which satisfies U(x) ~ U~oln(\x/a\), e.g. U(x) — 
0.5J7 ln[l + (x/a) 2 }. We will then use © to obtain the 
fluctuations of the time average showing that for high 
enough temperature the fluctuations increase with time. 
For this potential, for 1 < Ua/iksT) < 5, due to the slow 
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convergence of the tail of the distribution to P eq , and the 
slow power-law decay of E{x\,t) (which we shall shortly 
demonstrate), rendering £(xi) infinite for Uo/iksT) < 
2, one must consider the full time dependent problem 
instead of the time independent Eq. (|9]) and P eq (x). 
Generally the correlation function is given by 

/oo 
x 1 E(x u t 2 -t 1 )P(xi,t 1 \0,0)dx 1 . 
-OO 

(11) 

To solve this problem we used two approaches; the first 
is based on an eigenfunction expansion of the solution of 
the FP equation Such a calculation is lengthy and 
hence we adopt here a scaling approach. As seen from Eq. 
(TTT1) the key quantity to calculate is the ensemble mean 
E(xi , t) using Eq. ([8]) . Due to the homogenous character 
of the large x Fokker-Planck operator, it is natural to 
adopt a scaling ansatz: 



(12) 



where a and /3 are scaling exponents. Since for short time 
E(xi,t) ~ X\ we have g{y) ~ y for large y and a = /3. 
Inserting Eq. (|12l) in Eq. ([8j we find to leading order 
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where Uq = Uq/QcbT) is a key dimensionless parameter. 
To achieve a i-independent equation, we must have j3 = 
1/2, typical of Brownian motion. Then, 



g(y) 
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where M(a, b, x) [also denoted iFi(a;b;x)] is the Rum- 
mer M function [lij and we rejected a second so- 
lution in terms of the Kummer U function since it 
does not satisfy the boundary condition E(xi,t) — > 
when r — > oo (i.e., relaxation to equilibrium). The 
constant c\ is found by matching the solution in the 
y — > oo limit which corresponds to short times. Using 
M(a,b,x) ~ exp(a;)r(fe).T a_6 /r(a) and g(y) ~ y we find 
ci = {r(3/2)/r[(3 + (7o)/2]}(4L>)- c>0 / 2 . In particular, for 
long times, E{x\,t) ~ t~ u °/ 2 , so as we claimed, £,(xi) 
diverges for U <2. 

Steady state cannot be used to obtain the correlation 
function. To complete the calculation, we must have 
P(xi,ti\0, 0) which was recently obtained |l5[. The equi- 
librium PDF, since it decays as a power law P cq (x) oc 

|x| _c/o would give, for 1 < Uq < 3, (x(t 2 )x(ti)) = oo for 
t\ = t%. This is an unphysical behavior: at finite time 
one cannot have an infinite value for the correlation func- 
tion, since the particle cannot travel faster than diffusion 
permits. Specifically in the limit of long t\ we have [l5| 



P(xi,ti|0,0)~ P cq (£i) 
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FIG. 2: The aging correlation function Eq. (|16p perfectly 
matches numerical simulation of the Langevin Eq. (flj). 



Since T(a, 0) = T(a) as t — > oo, thermal equilibrium is 
reached. Nevertheless, for the calculation of correlation 
functions one must take into account the finite time cor- 
rection which is represented by the ratio of T functions. 

Aging correlation function. Inserting Eqs. (|14ll5p in 
Eq. dill) we find the non-stationary correlation function 
for the temperature range 1 < Uq < 3: 



(x(t 2 )x(h)) ~ (x 2 (h))f 0o 



where 
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The behavior in Eq. (|T6|) is very different than the sta- 
tionary case where the correlation function is a func- 
tion of the time difference t 2 — t\. In this tempera- 
ture regime the equilibrium mean square displacement 



diverges, (x 



' cq 



oo, while the time dependent solution 



Eq. (HU gives [3 (x 2 (h)} = a 2 c 2 {ADt 1 /a 2 )^- u '>)/ 2 , 
c 2 = 2(a/Z)lr(l/2 + U a /2)(3~U )}- 1 . We find /^(O) = 
1 which implies that C(ti,t±) — (x 2 (ti)) as it should. In 
the opposite limit t 2 3> t\, we obtain 



(.r(/ 2 ).r(/ J ))^r,(,-U S ))( | 



(18) 



withc 3 = (3/2-C/ o /2)0i : r(2 + [/ o /2)/3r(3/2 + C/o/2). In 
Fig. [5]we compare our analytical Eq. (|17[) with Langevin 
simulations showing excellent agreement for various mea- 
surement times. 

As mentioned we assume that the partition function 
function Z is finite and hence the steady state P° q (x) is 
normalizable. This excludes the well known Bessel pro- 
cess [l2| which can be mapped onto U(x) — Uo In |x| with 
its singularity at the origin. It is important to empha- 
size that (x(t 2 )x(ti)) ~ 1/Z depends on the shape of the 
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TABLE I: Scaling behavior of (x 2 (t)) and (x 2 {t)) for various 
values of fr = U /(k B T). 



potential in the whole space through Z . Hence for the cal- 
culation of the correlation function the regularity of the 
potential on the origin is vital. Interestingly this is not 
the case for all observables; e.g., E(xi,t) Eqs. (|12I14|) 
is Z independent and hence related to the Bessel process 
0. 

Ergodicity of the dynamics is classified in four domains 
which arc controlled by temperature, 
(a) The most interesting case is the regime 1 < Uq < 3. 
As we showed, a normalized steady state exists and from 
symmetry (x) cq — 0. If we naively assume ergodicity 
x — > (x) cq = and lim t _ J . 00 (:z; 2 (t)) = 0. Rather, from 
Eqs. ipTHjl we find 
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Changing variables to s = tz/tx — l^w = tjt\ — 1, we find 



(x 2 (t)) ~ c 4 (x 2 (t)) cx r 



(20) 



where c 4 =4 / °° diu(l + u;) (c>0 " 7)/2 (7- t/o)" 1 /^ M- We 
see that the fluctuations grow with time, hence ergodicity 
is broken. We find that c 4 w 0.2397 for fr = 1 and that 
it decreases monotonically to C4 = at Uq = 3. 

(b) For lower temperature, 3 < Uo < 5, the integrals 
in Eq. JTDJ still diverge, and (x 2 (i)) decays as t^-^o)/?. 
indicating an anomalously slow approach to ergodicity. 

(c) For Uq > 5 the temperature is low enough that Eq. 
([10)1 is now valid. For U(x) = 0.5U ln[l + (x/a) 2 ] we find 



(x 2 (t)) 



2(L>o - 4) 
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(21) 



which diverges when f7o — > 5. 

(d) Finally, for very high temperatures Uo < 1, the equi- 
librium state Eq. @ is not defined as the partition func- 
tion Z diverges. Here (x 2 (t)) cx t, exactly the diffusive 
behavior of a free particle, Uo = [l3j . 
These four different behaviors are confirmed via numer- 
ical simulations presented in Fig. [3J which illustrates 
convergence on reasonable computer time scales. A sum- 
mary of the scaling regimes is presented in Table 1. 

Relation with experiment. After the submission of this 
manuscript, an experiment on anomalous diffusion of 
ultra-cold atoms which employs the well known Sisyphus 
cooling scheme was reported [l(| ■ In the semi-classical 
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FIG. 3: In the non-ergodic phase k B T > U /3 (U < 3), 
{x 2 (t)) increases with time. For the critical point fcsT = Uo/S 
(Uo — 3) the fluctuations are constant. The dashed curves 
are theoretical predictions Eqs. (|20I21|I which agree very well 
with the numerical simulations. 



approximation, the atomic velocity distribution follows 
Fokker-Planck dynamics in an asymptotically logarith- 
mic potential [1, [TH, [I?} • Our work provides the theoreti- 
cal mean-square displacement in this experiment by iden- 
tifying our position x with the velocity v of the atoms. 
The measured atomic position is x(t) = Lv(t)dt and 
hence x(t)/t corresponds to the time averaged velocity. 
The PDF of the atoms in the experiment has been de- 
scribed with a Levy distribution, with a divergent vari- 
ance. However, our results show that the mean square 
displacement is finite for any finite measurement time. 
These seemingly contradicting findings are related to the 
well known dilemma whether Levy flights are at all phys- 
ical, since they predict divergin g m ean square displace- 
ment, which must be tamed 0, Ha]. We speculate that 
the Levy distribution found in the experiment describes 
the center part of the packet, which eventually is cut off 
to give a finite mean square displacement. Furthermore, 
using our results one can estimate the time in which the 
atoms remain within a finite domain, which is of course 
crucial for experiments. Experimentally one may also 
control the depth of the optical potential, here modeled 
with Uq and hence explore the nontrivial dependence of 
our results on this parameter. We will elaborate on these 
interesting points in a longer publication. 

Discussion. Aging correlation functions and ergodicity 
breaking typically describe glassy dynamics [ill [2(| (and 
Ref. therein). Our work shows that aging and ergodicity 
breaking can be found also for simple Markovian dynam- 
ics, without the need to introduce heavy-tailed waiting 
times into the kinetic scheme, nor disorder or many-body 
physics. The aging correlation function (|16p has a signa- 
ture very different than most previous work. The prefac- 
tor (x 2 (ti)) grows with time, and hence we call it super- 
aging. This is in contrast to normal aging where the cor- 
relation function is of the form C{ti,t\) = (x 2 ) qf(t2/ti) 
with a finite equilibrium value (x 2 ) eq . A similar non- 
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normal aging behavior, albeit with a logarithmic time 
dependence, has been found in Sinai's model of diffusion 
in a random environment fl9l ] . Unlike previous scenarios 
to ergodicity breaking, the amplitude of the stochastic 
process x(t) in our work increases with time, since the 
particle explores more and more of the tails of the equi- 
librium PDF as time goes on. Thus rare events where 
the amplitude x(t) of the Markovian process attains a 
large value are responsible for the non-ergodic behavior. 
This is clearly related to the power law tail of the equi- 



librium steady state P cq (x) oc |ir|~ c/o . More importantly, 
physical systems with fat tailed equilibrium states are 
common and hence this type of ergodicity breaking may 
find broad applications. 
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